Tunnelling-charging Hamiltonian of a Cooper pair pump at large Ej/Ec'- 
Modified Hamiltonians and renormalisability 
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The properties of the tunnelling-charging Hamiltonian of a Cooper pair pump are well understood 
in the regime of weak and intermediate Josephson coupling, i.e. Ej < Ec- Instead of perturba- 
tive treatment of charging effects, the present work applies the charge state representation in the 
the strong coupling case. From the discrete Hamiltonian we construct effective, truncated PDE 
Hamiltonians and analytically obtain approximate ground-state wave functions and eigenenergies. 
The validity of the expressions is confirmed by direct comparison against the results of numerical 
diagonalisation. For uniform arrays, our results converge rapidly and even (^-dependence of the wave 
function is described reasonably. In the inhomogeneous case we find the Hamiltonian to be para- 
metrically renormalisable. A method for finding inhomogeneous trial wave function is explained. 
The intertwined connection linking the pumped charge and the Berry's phase is explained, too. As 
addendum, we have explicitly validated the ground state ansatz for (f> — when N < 42. 



I. INTRODUCTION 

Josephson junction devices, e.g. Cooper pair boxes, 
superconducting single electron transistors (SSET) and 
Cooper pair pumps, have been extensively studied both 
theoreticallybViVi 6 and experimentallyi 7 i 8 i 9 i 10 i 1:L For a 
recent review, see Ref. 0. Possible applications in- 
clude at least direct Cooper pair pumpingjS, decoherence 
studies^ related metrological applications^ and the use 
of Cooper pair charge qubits or persistent-current qubits 
(SQubits) in quantum computationiiii 

The ideal tunnelling-charging Hamiltonian of a Cooper 
pair pump has been studied in detail in Refs. l2l ll llllfil 
Charge transfer due to direct supercurrent and adiabatic 
pumping due to varying gate voltages have been ade- 
quately described when the Josephson coupling is weak 
or at most comparable to charging effects. The case 
of strong Josehpson coupling in ideally biased arrays is 
still relatively unexplored. A single Josephson junction is 
known to be described by the Mathieu equation^ in the 
phase representation. For a superconducting single elec- 
tron transistor (SSET) the charge state representation is 
identical to one-dimensional discrete harmonic oscillator 
and, thus, the Mathieu equation in the island's phase 
representation>iZii& 

In this paper we first develop a method for obtaining 
an approximate solution of the Mathieu equation. Later 
on, we generalise the method for several dimensions and 
make the required corrections for our model Hamiltonian. 
In short, starting from the discrete Hamiltonian we con- 
struct a modified partial differential equation (PDE) for 
which a trial solution is obtained. Subsequently, the solu- 
tion is overlaid as the wave function the discrete Hamil- 
tonian and the result is compared against numerically 
obtained eigenstate. 

In order to sum up the obtained results we state the 
following: For homogeneous arrays of arbitrary length we 
find analytical and rapidly converging wave functions and 
eigenenergies. These expressions are derived from the 
developed method The case of non-zero phase difference 



is treated in a fairly satisfactory way. Inhomogeneous 
arrays are first treated by parametric renormalisation 
which yields an accurate approximation for the ground 
state energy. A modification of the original method im- 
proves the wave function, but not the asymptotical rate 
of convergence. 

Skeel and HardyiS have performed analysis on con- 
structing modified Hamiltonian when integrating systems 
of PDE's over time, see also Refs. I20II21II2I In these 
works numerical discretisation is approximately counter- 
acted by using a suitable truncation of the modified equa- 
tions. The principles of the present method are similar, 
although it is applied on a discrete eigenvalue problem. 

This paper is organised as follows. In Sec.[H]thc Hamil- 
tonian is defined and its structure is explained. In Sec. lllll 
we find an approximate solution for the Mathieu equation 
in charge state representation and postulate the general- 
isation of the method for several coordinates. In Sec. II VI 
homogeneous arrays are examined and explicit trial wave 
functions for the ground state are constructed. In Sec.lVl 
the developed formalism is extended to into account non- 
zero values of phase difference across the array. In Sec. lVIl 
the Hamiltonian is shown to be parametrically renormal- 
isable in the inhomogeneous case. Wave function is also 
constructed although the accuracy is not as good as in the 
homogeneous case. Finally, the conclusions are drawn in 
Sec. EH 



II. CONSTRUCTING THE HAMILTONIAN 

A schematic view of the system is shown in Fig.^ We 
assume that the gate voltages V g j are independent and 
externally operated. The bias voltage across the array, 
Vb, which controls the total phase difference 4> according 
to d<f>/dt = —2eV/h, is assumed to be ideally set to zero. 
Hence, <f> remains fixed and becomes a good quantum 
number in proper variables, which have been presented 
e.g. in Ref. |23|. On the other hand, a precise value of 4> 
means that its conjugate variable M, the average num- 
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ber of tunnelled Cooper pairs (M := —id/ dip), becomes 
completely undetermined. 

In the following, the tunnelling-charging Hamiltonian 

H = H C + H, ]l (1) 

is assumed to be the correct description of the micro- 
scopic system. The ideal model Hamiltonian simply 
neglects quasiparticle tunnelling as well as other de- 
grees of freedom. The most important parameters are 
the (average) Josephson coupling energy Ej and (av- 
erage) charging energy, defined as Ec ■= (2e) 2 /(2C). 
These determine "the Josephson-charging ratio" which 
is denoted by £j := Ej/Eq- The Hamiltonian and 
the operation of a Cooper pair pump in the weak cou- 
pling regime is further determined by (normalised) gate 
charges q := {q x , . . -,q N -i}, where q k := -C^ k V^ k /2e. 
In the present model relative junction capacitances c, 
where c k := C k /C and J2k=i ^k, 1 = N/C, also deter- 
mine individual Josephson energies by Ej^ k '■= c k Ej. For 
uniform or homogeneous arrays we have c k := 1, while 
the inhomogeneity can be reliably quantified by the in- 
homogeneity index X inh := [J2 k ( c k 1 - l) 2 /N}^ 2 M^ 




FIG. 1: An ideal superconducting array of independent 
Josephson junctions. Here C k and Ej $k are the capacitance 
and the Josephson energy of the k th junction, respectively. 
The total phase difference across the array, <j>, is a constant 
of motion. 

The matrix elements of the charging Hamiltonian Hq 
are given by the capacitive charging energy and thus they 
readi^ 

21 



(n\H c {q)\n)$ = E, 



c 



N , 



fc=i 



Ck 



/ N 
1 / Ufe 



N 



E 

\k=l 



Ck 



(2) 



where the number of Cooper pairs on each island is given 
by n 



(ni, 



,njv-i)- The quantities {vk} k=1 are an 
arbitrary solution of the charge conserving equations 

Vk - v k +i = n k - q k . (3) 

Tunnelling of one Cooper pair through the fcth junction 
changes \n) by 5 k , where the non-zero components are (if 
applicable) (S k ) k = 1 and (5 k )k-i — ~ 1- The tunnelling 
Hamiltonian is given by 

N 



Hi 



E 

n.k— 1 



CkEi 



(\n + 6 k )(n\e i 't>/ N +K.c.). (4) 



The supercurrent flowing through the array is determined 
by the supercurrent operator 7g = (—2e/K)(dH/d<f)), a 
Gateaux derivative 2 ^ of the full Hamiltonian. By chang- 
ing the gate voltages adiabatically along a closed path 



r, a charge transfer Q t 



Q p is induced. The 



pumped charge, Q p , depends only on the chosen path, 
while the charge due to direct supercurrent, Q s , also de- 
pends on how the gate voltages are operated. If the sys- 
tem remains in a adiabatically evolving state |m), the 
total transferred charge, Qtot> in units of — 2e, reads 2 *^ 



drjm(t) 



2 <j> Re 



(m\M\dm) 



(5) 



where \dm) is the change in \m) due to a differ- 
ential change of the gate voltages dq and rj m = 
— Jq (E m (t)/H)dt is the dynamical phase of the wave 
function. 

Clearly, the pumped charge is closely related to the 
the geometrical Berry's phased j m (T) — i <f> r (m\dm) . 
The pumped charge can be evaluated from Eq. (J5J in the 
charge state representation once the overall phase of the 
eigenstate is fixed consistently for all q. If the examined 
state is sufficiently non-degenerate for all values of 0, the 
eigenstate can be expanded as a Fourier series in with 
real coefficients {a^ ;}. Consequently, for a fixed value 
of (j> the differential pumped charge is given by a gauge- 
invariant expressio: 



,15 



dQ p (<p)=J2 



I'—O n.l — — oo 



2(1 + Y n /N) 



d(ctfi,i a ft,i+i') 



1 + <Vo 

+ l\anjda ri j + i, -an,i+i<dan,i)}cos{l' (f), (6) 



where Yfi is an additional class label. In constrast, a 
differential change in the phase difference <fi for fixed gate 
charges q induces no pumped charge, because we find 



dQ p (dcj>) = 21m (Mm\Mm) 



d<p = 0. 



(7) 



Now consider the Berry's phase 7 m induced by an in- 
finitesimal closed cycle C at (q, (j>) with sides dq and d<f> 
as shown by the l.h.s. pf Fig.|2 The result divided by d<f), 
i.e. djm /dcj), is identical to dQ p apart from the sign of 
the first term. In other words, the contribution from the 
first and third part of the cycle gives the non-integrable 
part of dQ p , while the second and fourth part add up to 
the integrable part multiplied by —1. Thus the path for 
which the (^-"derivative" of Berry's phase is identical to 
dQp{4>) is n °t a closed cycle but a more complex path 
illustrated in the r.h.s. of Fig. [3 

From here on the expression for the charging energy, 
Eq. (J2J, is examined in detail. This is done in order 
to rewrite the Hamiltonian in as simple a form as pos- 
sible. In the homogeneous case the quadratic form is 
easily diagonalised and we find N — 1 identical eigenval- 
ues of Ec and one zero-energy mode in the direction of 
vq := (1,1,..., l)/\/~N. This demonstrates the unique- 
ness of the charging energy expression for each charge 
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FIG. 2: An infinitesimal cycle C corresponding to Berry's 
phase 7 m (C) consists of four legs. The charge transfer Q p for 
a fixed (j> is identical to the Berry's phase induced by travers- 
ing the legs in the shown directions. This path can not be 
followed continuously in the (q, 0)-plane. 



state and, consequently, the same zero-energy mode is 
observed in the inhomogcneous case, too. 

In a proper representation of the g-space, the charg- 
ing energy for homogeneous arrays can be expressed as 
-^cllsllii where || • H2 is the usual Euclidean norm. Thus, 
the representatives of the tunnelling vectors {Sj}, de- 
noted by qj, are required. Above all, they must be nor- 
malised according to 



qj ■ q k = Sj k - l/N. 



(8) 



In an orthonormal (N — l)-dimcnsional basis, where ej ■ 
efe = Sjk and x — (xi,X2, . . . ,xat_i), the representatives 
define variables {qj} according to 



JV-l 



Qj 



(x) := ^ ( e k ■ qj)xk- 



(9) 



k=l 



The normalisation condition JSJl yields relations 



JV 



JV 



JV-1 



(10) 



3=1 



3=1 



3 = 1 



which are valid for all values of x. 

Suitable representatives for cases N = 3 and TV = 4 
are easy to find and their visualisation is obvious. When 
N = 3, we select 

qx = (v/273,0), & = (-I/a/6, 1/V2), 

q 3 = (-1/V6,-1/V2), (11) 

which describes three directions separated by identical 
120° angles. The resulting transformation of coordinates 
and the so-called honeycomb structure is shown in Fig. [31 
The gate charges q\ and c?2 determine the origin of the 
induced, rectangular coordinate system (x\,X2)- 

For N = 4, symmetric representatives are given by the 
well-known body centered cubic lattice (BCC) of solid 
state physics, explicitly 



(1,1, -1)A 
(-l,l,l)/2, 



<74 



(1,-1, l)/2, 
(-1,-1, 



(12) 




FIG. 3: On the left-hand-side, the so-called honey comb 
structure induced by the charge state lattice (N = 3). The 
regular lattice determined by the representatives q\, 52, and 
<?3 is shown on the right-hand side. The origin of a new rect- 
angular coordinate system (xi,X2) is set by the gate charges 
q — (31,92). The charging energy for a charge state (gray 
circle) then reads Ec(xi + X2), where the nearest-neighbour 
distance has been scaled to y2/3. 



These representatives are convenient when studying the 
case = 4, but a more general method for obtaining 
representaves is required. By augmenting the existing 
representatives for A^ we can always obtain the set for 
N + 1. The additional representative is set to lie along 
the new (first) coordinate axis with the correct length 
\/N/ (N + 1). The normalisation condition (JjJJ is satis- 
fied if all other representatives are retained as they were 
with an identical first component of —l/y/N(N + 1). 

Applying this method inductively, starting from trivial 
case of A" = 1, yields the general representives for any 
N. Let the length of the array be N and denote the 
jth re p reseir t;ive and its fc th component by q^ and , 
respectively. The A^ — 1 components are obtained from 
three simple rules: 



n N - 
n N 



y/(N-j)/(N+l-j) 



(0 

(«) 

(»*) 4k<3) = -l/y/(N-j)(N+l-j) (13) 

The above transformation simplifies and symmetrices the 
tunnclling-charging Hamiltonian for arrays of any length. 
Inhomogeneous arrays can also be considered once the 
tools have been developed. 



III. MATHIEU EQUATION AND DISCRETE 
HARMONIC OSCILLATOR 



The canonical form of the Mathieu equation reads^ 

(14) 



dv 



^ + {a-2qcos(2v))y = 0, 



where y(v) is the solution, q is a parameter and a is 
known as the characteristic value or eigenvalue. 

The Hamiltonian of a SSET for a fixed phase difference 
cj> in the charge state representation can be mapped onto 
a one dimensional discrete harmonic oscillator (DHO), 



4 



see e.g. Our chosen form includes a nearest- 

neighbour coupling — ej/2 and the potential V(n) — 
(n — no) 2 , where n is an integer. The equation for the m 
amplitude a n now reads 

(n - n ) 2 a n - (ej/2)(o„_i + a n+1 ) = Ea ni (15) 



where E is the eigenvalue we are looking for. In order to 
obtain the solution of the discrete equation, we assume 
that a n is a continuous function an replace other ampli- 
tudes by respective Taylor expansions. We denote the 
step size by h (here h = 1) which yields 



-h + O-n+h = 



2k j2k. 



2h 2k d 



k=Q 



(2k)\ dn 



2k 



(16) 



a differential equation for a n . 

We now transform into conjugate variables of the is- 
land charge, i.e. n — no — > —id/dO and id/dn — > 9. Col- 
lecting the terms, we find 



d 2 a(9) 
d9 2 



[E+(e 3 /2)cos(9)}a(9)=0, 



(17) 



which is identical to the Mathieu equation with a = 4E 
and q = 2ej once we choose 9 = (9 + 7r)/2. 

In the limit q — > oo, the ground state energy can be 
read from Eq. 20.2.30 of Ref.El with the result 



E, 



DHO 




V2" 



1 

16" 



VV^i 3/ej 
256 2048 



+0(ej 3/2 ), (18) 



confirmed by numerical diagonalisation, too. Returning 
to the charge state representation, we divide the eigen- 
value problem by ej and define the oscillator frequency 

The lowest 



lo := y/2/ej and scaled energy E := E/ej 
order approximation becomes 



■1. 



-a" /2 + \u?{n - n„) 2 a = E, 



(19) 



which is analytically solvable with E = lu/2 and a(n) oc 
exp(— u>(n — no) 2 /2). From here on, no is omitted from 
the expression [n — no) for brevity. 

The discretisation naturally affects the wave function 
and as well as the eigenenergy (|18|) . The lowest order 
approximation ipi for the discrete wave function is 
naturally a Gaussian wave function. The optimal, but 
unnormalised, wave function is given by 



i\)\{n) oc exp 



un 2 /2 
1 - w/8 



tpi x (n) oc n exp 



un 2 /2 
1 - 27w/16 

More accurate a wave function reads 



, wn 2 (l-w 2 n 2 /48)/2\ 
Woe exp ( j 



(20) 



(21) 



(22) 



n) oc 



24 



exp 



o;(n 2 -cj 2 n 4 /48)/2 \ 
1 - 147^/640 J 



where a cutoff (1 — ui 2 n 2 /48) — > 1 must be applied for 
large enough values of n, when the deviation becomes 
greater than 20-30 %. These wave functions have been 
compared against the result of numerical diagonalisation, 
ij)d, by taking the norm of the difference, in short \\ipd — 
ipj\\, which yields approximately 



0.018/v^J, 
0.009/ej. 



(24) 
(25) 



Because both trial wave functions converge towards the 
actual eigenstate of the system, approximate eigenener- 
gies corresponding to ip\ and ip2 can be easily evaluated. 
Setting no = and examining the equation for coefficient 
a gives 



-£j(oi/oo), 



(26) 



where £j — > oo. Expanding the terms in powers of lu = 
\/2/ej gives the desired result. We find that E^ !1 first 
deviates from the constant order in which the term is 
— 1/8 instead of the correct —1/16. As expected E^, 2 is 
much better, and even the term — a/2/ej/256 is correctly 
reproduced. 

The significance of the corrections in ^{n) with re- 
spect to the continuous solutions is relatively clear. The 
denominator 1 — 3w/16 cancels 3/4 of the of the lead- 
ing second order term — (cj/2) 2 /2 and gives the correct 
eigenenergy in the constant order. On the other hand, 
the term proportional to n 4 is related to the truncated 
differential operator 



dn 2 



1 d 4 
Vldr?' 



(27) 



The coefficient -lj 2 /48 can be divided in two parts, 
namely 1/12 and — cj 2 /4 which seems reasonable as the 
latter scales correctly as function of oj, while the former 
changes if the step length h is altered. 

This approach is rather similar to that of Skeel and 
Hardy 19 although they consider time-dependent prob- 
lems instead of eigenvalue problems. Systems of dif- 
ferential equations are replaced by modified equations 
which try to compensate for the discretisation error. The 
present potential is harmonic and in the conjugate rep- 
resentation truncated potentials are anharmonic in na- 
ture. Anharmonic oscillators have been studied, and ex- 
act eigenvalues have been obtainedi^L 2 ^ Unfortunately, 
the sign of our leading correction is negative, so these 
works are not applicable here. 

If there are two orthogonal and independent directions, 
the wave function factorises and the one-dimensional re- 
sult can be generalised. Nevertheless, we make the fol- 
lowing assumption which is to be justified later. Let 
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our Hamiltonian be defined on a regular, discrete lattice 
of the coordinates x and the potential be isotropic and 



harmonic, i.e. V(x) 



/2. Interactions between 



(neigbouring) lattice sites are expanded in terms of par- 
tial derivates up the fourth order in a similar manner to 
Eq. Ijl6|l . We postulate an analytical trial solution if the 
second order operator is the Laplacian and the modified 
PDE eigenvalue problem has the form 

-i (V 2 + D 4 /12) ip + ±uj 2 \\x\\ 2 iP = EtJj, (28) 



where D4 is a fourth order partial differential operator. 
We define corresponding 'conjugate variable' D4 by re- 
placing each partial derivative with respect to Xj by Xj 
itself. For example, if D4 = (V 2 ) 2 , i.e. the square of the 
Laplacian, the conjugate variable is the fourth power of 
the norm, explicitly, D4 = (||i?|| 2 ) 2 = \\x\\ . Our unnor- 
malised trial wave function is given by 



ip2{x) oc exp 



-(w/2)(||x|| 2 -w 2 £> 4 /48) 



1 



(29) 



where a is chosen so that cancellation of the constant 
order term in energy is exactly 3/4, just as the factor 3/16 
in Eq. (|22|l . The conjugate variable D4 gives the correct 
functional form, although a cutoff for too large values as 
compared to \\x\\ 2 must be naturally applied. Hopefully, 
the asymptotic convergence of the norm — 1P2W is 
better than l/^fej. The general asymptotical solution 
for the discretised harmonic oscillator has geen recently 
given in Ref . I29T 



the eigenvalue problem into 



1 N 

-Y 

2 ^ 

3=1 



{ag-Sj- 2a 3 +a s+qj ) + ^u 2 \\x\\ 2 a s = Ea Sl (31) 



Using the procedure explained in the previous section, 
we can find the corresponding modified PDE. The trun- 
cation means that each term (ag-qj — 2aj + ag+q j ) cor- 
responds to a second order derivative and a fourth order 
derivative. The sum of the second order derivatives yields 
the Laplacian operator V 2 due to the second part of 
Eq. ijTUI) , and the form of the modified equation matches 
Eq. (|28p. Next, we must evaluate the form of D4, find cor- 
rect value of a, and compare the resulting wave function 
and eigenenergie against numerically obtained results. 
The simpler, optimal Gaussian wave function reads 



■01 {x) oc exp — 



2 /2 



1 - (JV- l)w/(8AT) 



(32) 



where the modification of (N — 1)/N in the denominator 
arises from the fact that ||^|| 2 = (N - 1)/N. The op- 
timally as well as the expected rate of convergence, i.e. 
has been confirmed up to N = 10. 
In case N = 3 we find that 



D 



(JV=3) 



(V 2 ) 2 /2 



(33) 



(N=3) 



and hence D 
improved wave function for N = 3 reads 



72. Because ||gj|| 2 = 2/3, the 



IV. STRONG JOSEPHSON COUPLING AND 
HOMOGENEOUS ARRAYS AT <j> = 

When the Josephson energy Ej is large as compared 
to charging energy Ec it seems preferable to express the 



Hamiltonian in terms of the phase differences 



We 



choose to remain in the charge state representation for 
two reasons. First, the charging Hamiltonian is diffi- 
cult to evaluate in the independent phase representation. 
Additionally, the model Hamiltonian is already diagonal 
with respect to the total phase difference </>. The model 
Hamiltonian can be approximately diagonalised and in- 
teractions between states with different values of <fi should 
be included later. 

The Hamiltonian equation is explicitly written in units 
of Eq, and more specifically, each equation (row) of the 
eigenvalue problem is examined separately. Each charge 
state n is labeled according its position in the orthonor- 
mal coordinates x — (x\, . . . , xn-i)- In units of Eq, the 
equation for the coefficient as reads 



^lE^^+^-ft) = Eas - (30) 



TV 



Now, consider the case = and large values of ej in de- 
tail. Writing the eigenvalue as E := N + E/ej transforms 



/ (^=3) 
lp2 ex P 



HHI72)(i-^ 2 IHI796) 

1 - u/8 



(34) 



This proves to be quite accurate as the norm of the error 

vanishes according to \\ipj — i/^ -3 ^ I ~ 0.0045/ej. As 
arrays become longer, pure radial (energy) dependence 
is not enough, since the operators D4 become more com- 
plicated. For the BCC representatives (N = 4) the dif- 
ferential operator is given by 



Da 




1 3 



o 4 



3=1 J 



corresponding to a wave function ijj^ proprotional to 



(35) 



/ 



exp 



(w/2) 



INT- V 3 -i 

64 i-ii=\ 96 



1 - 9w/64 



(36) 



Because (x 4 ) — (||x|| 4 )/5, this also explains why the best 
energy dependent fit occurs at 3w||x|| 4 /320. 

For longer arrays the expression for the fourth order 
differential operator becomes quite complicated and less 
informative. Fortunately, the value of the conjugate vari- 
able D4 can be easily obtained for any point x. The 
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simple expression is based on inner product of the x-space 
and the representatives qj , in short 



A' 



(37) 



The differential operator Dn can be read from the above 
expression by retaining the components of x in symbolic 
form and transforming each coordinate its corresponding 
partial derivative. The correct cancellation requirement 
implies that the general form of a is given by a = 3(N — 
1)/16TV. 

Thus, the general trial wave function ij)^ becomes 



-02^ (%) = A exp 



(u;/2)(\\x\\ 



co 2 D { 4 N) (x)/48) 



1 - 3(N - l)w/(16iV) 



(38) 

where A is a normalisation factor and Dn{x) given in 
Eq. I|37|) is evaluated for all charge states in the used ba- 
sis. A suitable cutoff with respect to uj 2 Dn(x)/(A8\\x\\ 2 ), 
e.g. between 0.2 and 0.3, is naturally important. The 
wave function is independent of the representives {qj}, 
but those given in Eq. (fH!)l are probably the most conve- 
nient. The rate of convergence of the norm \\ipd — V^'ll 
has been confirmed as 1/ej up to TV — 7. Tentatively, 
the same applies for N = 10, although diagonalisation 
was limited below £j « 20. 

The ground state energy is virtually independent of 
the gate charges q when £j is large enough. Thus 
can be approximately obtained as in the one-dimensional 
case, see Eq. <|26ll . All 2N neighbouring amplitudes are 
identical which now gives 



(AT) 



-Nej exp 



( W /2) 



N-l 
N 



1 - 3(TV-l)w/(16TV) 



(39) 



Expansion in powers of u yields the asymptotic expansion 



EW~-Nej + (N-l) 



16N 



O 



(40) 

verified by direct comparison against the numerically ob- 
tained eigenvalue for cases which allow diagonalisation. 
No analytical expression for the term proportional to 
have been found, but it is not correctly repro- 
duced, either. Direct calculation, using a method pro- 
posed in Refs. and |3(], validates the above ansatz 
and corresponding asymptotical eigenenergy for TV < 42, 
thought We now proceed to the the case when is no 
longer zero. 



neighbour coupling contains a term e ±lc ^/ Ar . When <ft is 
sufficiently small the phase does not vary significantly be- 
tween nearest neighbours and as the first approximation 
the phase can be neglected in the corresponding equa- 
tions. We then consider the absolute value of the ampli- 
tudes and observe that the differential operator is simply 
multiplied by a factor cos(</>/TV). 

Consequently, the approximate eigenvalue problem to 
the original one, except that lu is replaced by w — 
u> I d cos(</> / TV) . The ground state energy can be obtained 
from Eq. (|40|l with £j — > £j cos(0/TV). The accuracy of 
this expression is rather good, even for large values of <j) 
if £j is sufficiently large. The convergence in terms of 
the absolute values of the amplitudes is satisfactory, too. 
Convergence in terms of trial wave function goes 
clearly as and that of l^l goes nearly as l/e,j, 

weakening as 4> increases. 

In order to consider the complex wave function explic- 
itly, the approximate differential operator induced by <p 
must be constructed. The first order differential opera- 
tor is always cancelled on behalf of the first property in 
Eq. l|10fl . The common prefactor of the third order terms, 
relative to the Laplace operator, is here — i sin(0/TV)/3. 

Because the conjugate coordinate D± ' was so success- 
ful in describing the homogeneous case, we define a third 
order conjugate coordinate which evaluates to 



N 



(41) 



The first guess for the phase of the trial wave function is 
then given by 

M*) „ f (u>/2)(-iu>D^)/(6N) \ 

where sin(0/TV) has been approximated by <p/N . Numer- 
ical diagonalisation clearly confirms the dependence on 

~ (3) 

D N , although a numerical correction factor of the or- 
der of 0.7-0.75 for all TV has to be added. Additionally, 
but expectedly, the phase dependence is slowly damp- 
ened for larger values of ||a?||. Yhe magnitude of these 
amplitudes rapidly decreases which makes the imaginary 
components even smaller. Thus, the leading component 
of the phase simplifies to 



exp(i& sin(0/TV)w 2 D^ ) / 12 ): 



(43) 



where b, 



4> 



0.7. Finally, we turn in the direction of 



inhomogeneous array. 



V. EFFECTS DUE TO NON-ZERO PHASE 
DIFFERENCE 

For non-zero values of the phase difference <p the wave 
function becomes complex valued because the nearest 



VI. INHOMOGENEOUS ARRAYS AND 
RENORMALIS ABILITY 

Our main aim is to obtain a wave function similar to ip\ 
in the inhomogeneous case at (f> = and, subsequently, 



7 



improve this wave function. Effects due to non-zero 4> are 
treatable in principle, but the expression become rather 
messy and accuracy is not that good, ft suffices to say 
that the behaviour of the eigenenergy corresponds to the 
effective coupling strength ej cos(4>/N). 

In the inhomogeneous case the charging energy reads 



Eq 



N „2 1 / N 



(44) 



The biasing to zero voltage implies that ^2jLi ( v j/ c j) = 
0, although the above expression is invariant under trans- 



formation Vj 



■ y. As each coupling is multiplied by 



the second order approximation for the Hamiltonian 



3 • 

becomes 



N 



+ 



2 



2^/32^)2 2c 



(45) 



the homogeneous value is dominantly proportional to the 
square of the inhomogeneity index X^. The result, 



-(N If (2N 3)Xf nh + £? =1 (&f/iV) 



16N 



(50) 



has been confirmed up to N — 6 if only a single ca- 
pacitance deviates from the others. In case N = 3 this 
expression has been tested more rigorously and further 
corrections do vanish as 1/ y/gj. 

In order to improve the results, more elaborate trans- 
formations are required. The most viable transforma- 
tion is based on diagonalising the charging energy and 
transforming the representation space (x-space) in such 
a manner that the charging energy is proportional to the 
square of the new norm. New representatives q'j are ob- 
tained and the differential operators in the second and 
fourth order can be obtained. For some special cases, 
the second order differential operator is of the Laplace 
type, i.e. the conjugate coordinate is given by 



where j3 = \f N/ (N — 1). For sufficiently small values 
of u> and reasonably homogeneous arrays the condition 
^2jLi(vj /cj) = does not vary much between neighbour- 
ing points. In other words, the error between different 
lines of the eigenvalue equation is insignificant. Under 
those circumstances we renormalise the coordinates ac- 
cording to 



N 



J • 



(46) 



which yields a Hamiltonian identical to the homogeneous 
case. In a similar manner, we write the the lowest order 
wave function as 



1-(N- 1W8 



ip^ nhS> (x) cx exp ( — 



(47) 



where the summation gives simply the charging energy 
corresponding to x. This is the best Gaussian wave func- 
tion in the renormalised coordinates w, and the rate of 



convergence of the error the expected 1 / y/ej . 

The Hamiltonian of an inhomogeneous Cooper pair 
pump is thus renormalisable and the leading terms in 
the eigenenergy are 



E 



(AT) 
0,inh 



N 



(N 




Oil). 



(48) 



The constant term can also be evaluated if we assume 
a cancellation of 3/4 in this term which is correct for 
homogeneous arrays. We simplify the expression 



1 N 1 



16 ^— ' Cj 

.7 = 1 J 



(49) 



by denoting bj := (1/cj — 1) and collecting the terms. 
Not so unexpectedly, and as in Ref.Q the deviation from 



j'=i 



(51) 



In those cases the fourth order coordinate 

N 



i=i 



(52) 



yields a trial wave function which can be compared 
against the numerically obtained wave function. In most 
cases, the Laplacian operator is slightly distorted, but 
for small inhomogencities this can be neglected as the 
first approximation. In both cases the results are not as 
good as in the homogeneous case, but the improvement 
with respect to Eq. (|47fl is significant. Due to dimensional 
limitations the comparisons between wave functions have 
been performed when N = 3. 

As shown by the cancellation in the eigenenergy, no 
isotropic value of a such as 3(A^ — l)/16iV in Eq. 138f) is 
can be used. Rescaling of the coordinates changes the op- 
timal value of a in different directions, and some further 
improvement may be obtained by using a non-isotropic 
a(x) in the calculations. Minor improvements can be ob- 
tained by fiddling with the coefficients of the coordinates, 
too. We conclude this section by stating that significant 
improvement of the wave function has been obtained, but 
so far no analytical expressions have been able to reach 
asymptotical convergence better than 1 / ^/ej. 



VII. CONCLUSIONS 

We have developed a method for obtaining an (approx- 
imate) analytical solution for Laplace type eigenvalue 
equations with a harmonic potential and discreteness in- 
duced higher order corrections. In the one-dimensional 
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case corresponding to the Mathieu equation the re- 
sults were convincing and thus we applied the proposed 
method on the tunnelling-charging Hamiltonian of an 
ideally biased Cooper pair pump. 

We have obtained reliable analytical expressions for the 
ground state wave function and energy for homogeneous 
arrays of arbitrary length. Furthermore, effects due to 
nonvanishing phase difference were relatively well de- 
scribed and the Hamiltonian of an inhomogeneous pump 
was shown to be renormalisable. Again, reliable eigenen- 
ergies and reasonable eigenfunctions were obtained. Fur- 
ther improvements in the inhomogeneous case have been 



proposed and partially carried out, too. 
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